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Abstract: The analysis and proper documentation of the properties of 
closed-loop control software presents many distinct aspects from the analysis of 
the same software running open-loop. Issues of physical system representations 
arise, and it is desired that such representations remain independent from the 
representations of the control program. For that purpose, a concurrent program 
representation of the plant and the control processes is proposed, although the 
closed-loop system is sufficiently serialized to enable a sequential analysis. While 
^ dealing with closed- loop system properties, it is also shown by means of examples 

how special treatment of nonlinearities extends from the analysis of control 
specifications to code analysis. 



^ 1 Closed-loop software system analysis 

00 

The first part of this report essentially dealt with the analysis of open-loop 
control systems, that is, the analysis of control algorithms without regards for 
^SJ their stabilizing effect on the plant: In |FA08| . we were concerned with verifying 

T-H and quantifying the boundedness of all computed quantities in the controller, 

00 regardless of the actual performance of the closed-loop system. In that regard, 

the analysis is similar to that described in [CCF+OS , although our work is much 



^> more focused on control algorithms and, unlike the work in |CCF+05] , it ignores 

the many other software components present in a fully developed avionics sys- 
'V^ tem. In the present document, we focus on the analysis and verification of the 

system's closed-loop performance, beginning with closed-loop system stabiHty. 
One of the elements of this report is the modeling, by means of an example, 
of closed-loop control systems as two concurrent processes, where one process 
represents the real-time computer program and the other represents the sys- 
tem being controlled. The concurrent representation of these processes does 
not make their analysis more complicated, and it provides the engineer with 
greater flexibility by making the choice of controller representations and imple- 
mentations independent from the representation of the plant being controlled. 
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Again, the best way to perform the closed-loop system analysis is to begin with 
specification-level closed-loop system analysis, and then to move on to incorpo- 
rate more detailed controller implementations. This is the process illustrated in 
this report. 



1.1 Specification-level analyses: Stability 



The analysis of the specifications of the lead-lag controller designed and pre- 
sented in |FA08| begins with Figure [l] In that figure, the controlled system 
is represented by a continuous differential equation, while the controller is a 
discrete-time system that interacts with the continuous system by means of a 
signal sampler at the input and a zero-order hold at the output. The dynamics 
of the controller may be written 

Xc,k+i = AcXc,k + Be SAT(?/fc) 
Uk = CcXc,k + Dc SAT(j/fc) 



with 



Ac 



0.4999 
0.0100 



-0.0500 
1.0000 



[564.48 0] , and 



Bc^ 
-1280. 



It is well-known from control theory |FP80j that the analysis of such a closed- 
loop, sampled-data system can be performed by replacing the plant, the zero- 
order hold and the sampler by the discrete-time system 



Xp,k+l 

Vk 



ApXpk + BpUk 

CpXpk 



with 



1.0000 0.0100 
-0.0100 1.0000 



1 Br, 



0.00005 
0.01 



and Cp = [1 0]. Moreover, Xp^k — a;p(0.01fc), where fc S N. Then the analysis 
of the system shown in Fig. [ijis equivalent to the analysis of the discrete-time 
system shown in Fig. [2j The stability analysis of this discrete-time system 
is made more difficult than that of the controller alone |FA08j because of the 
presence of a saturation nonlinearity aimed at limiting the range of sensor inputs 
to the control system. This nonlinearity may be accounted for in the overall 
system stability analysis by using a variety of computational techniques, such 
as those described in [ BEFB94 . In particular one of the available techniques 
consists of computing a quadratic Lyapunov function for the system. It is, for 
example, possible to show that the quadratic function 



V{x) 



Xii 



p 



Xc 
Xxi 



, with P 



0.2205 0.0188 -0.0750 0.0177 

0.0188 0.4736 0.0535 0.0015 

-0.0750 0.0535 0.1012 -0.0049 

0.0177 0.0015 -0.0049 0.0015 
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Xc,k+1 — 
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Uk = 


564.48 0] Xc,k - 1280 SAT(yit) 
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Figure 1: Feedback system 



^p,k-\-l — 


1.0000 
-0.0100 
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■ 0.00005 ■ 
0.01 


Uk, Xo = 


Xo 
Xo 




Vk = 


1 0]xp,fc. 
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Figure 2: Equivalent discrete system 
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decays along all trajectories {xc,k, Xp^k)k=i,2,... for initial conditions where the 
controller is at rest (xcfi = 0) and the initial state Xp^ of the mechanical system 
lies inside the ellipse £q with 



0.1012 
-0.0049 



-0.0049 
0.0015 



and £q = {y e R2 I y^Qy < 1} 



(1) 



The mechanisms used to reach that conclusion rely on standard stability con- 
siderations, most notably absolute stability theory: Consider the closed-loop 
system dynamics 

Xk+i ^ Axk + BSAT{Cxk), 

with X — \xT x'l]'^ , 



A = 



Ac 
BpCc 





An 



B = 



0.4990 -0.0500 

0.0100 1.0000 

0.0282 1.0000 0.0100 
5.6448 



-0.0100 1.0000 



Be 

BpDc 



and C= [0 Cp]. 



Consider the set of x's such that x^ Px < 1. It is easy to show that (i) the set of 
all X = [0 x'^y such that Xp e £q is contained in £p, and (ii) that for all x E £p, 
(SAT(C2:) - 0.2C2;)(SAT(Cj:) - Cx) < 0. Introducing yc,fc = SAT(Ca;fc), it is 
therefore sufficient to show that 



X h 



whenever {yc,k — 0.2Cxk){yc,k — Cxk) < 0. An equivalent condition is 

T 



xlPxk < 



X 
Vc 



B^ 



P[A B] 



P[I 0] 



X 
Vc 



< 



whenever {yc — 0-2Cx){yc — Cx) < 0. Using the well-known 5-procedure |AG64] . 
a sufficient condition is the existence of A > such that 



A^ 
B^ 



P[A B] - 



F [/ 0] - A 



0.2C'^C 
-0.6C 



-0.6C^ 
1 



< 0. 



(2) 



The latter inequality can be easily shown to be true for A = 6.76. It must 
be kept in mind that all statements related to numerical objects and compu- 
tations are made under the assumption that computations on reals are exact. 
This assumption is obviously false and must be eventually addressed when im- 
plementing the operational software verification tools that may arise from this 
paper, using for example techniques developped by Goubault [GouOlj . We are 
now interested in showing how this proof of (local) closed-loop stability may 
be translated at the code and system level. To be more precise, we will show 
how invariancc of the ellipsoid £p can be exploited to develop a proof of proper 
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"/o Controller Dynamics 



Ic: 
2c: 
3c: 
4c: 
5c: 
6c: 
7c: 
8c: 
9c: 
10c 
11c 
12c 



Ac = [0.4990, -0.0500; 0.0100,1.0000] 

Cc= [564.48, 0] ; 

Bc=[l;0] ; D = -1280; 

xc = zeros (2,1) ; 

receive (y) ; 

while (1) 

yc = max(min(y,l) ,-1) ; 

u = Cc*xc + Dc*yc; 

xc = Ac*xc + Bc*yc; 

send(u) ; 

receive (y) ; 

end 



% Plant Dynamics 



2p 
3p 
4p 
5p 
6p 

7p 

Bp 
9p 



Ap = [1. 0000, 0.0100;-0. 0100,1. 0000] ; 
Cp=[l,0] ; 

Bp = [0.00005; 0.01] ; 
while (1) 

y = Cp*xp; 

send(y) ; 

receive (u) ; 

xp = Ap*xp + Bp*u; 
end; 



Table 1: Concurrent representation of plant and controller implementations 



behavior, such as stability and performance, for the computer program that 
implements the controller, as it interacts with the physical system. Unlike the 
developments relative to controller open-loop properties, this proof necessarily 
involves the presence of the controlled artifact. It would not make sense for 
the system to be represented at various degrees of computer program imple- 
mentations since these do not have any physical meaning. Thus we choose to 
represent the plant and the computer program as two concurrent programs, as 
shown in Table [T] The computer program representation of the mechanical sys- 
tem will remain invariant, written in a high-level language like MATLAB, while 
the representation of the controller will be allowed to evolve to reflect the various 
stages of its implementation, without changing the representation of the phys- 
ical system. We assume the two programs communicate by means of send / 
receive commands. It is assumed that the receive command is blocking, that 
is, the program execution cannot proceed until the variable of interest has been 
received. While the programs in Table [T] are purely event-driven, it is possible 
to modify the algorithmic description of the plant to account for the presence of 
time. The question of establishing proofs of stability of the closed-loop system at 
the code level is necessarily tied to understanding the behavior of the controller 
and that of the plant. The entire state-space then consists of the direct sum of 
the controller's state-space and that of the plant. We now propose to leverage 
the approach developed in the first part of this report |FA08j to document the 
corresponding system of two concurrent programs. This documentation can be 
extended to other representations of the control program implementation: Ta- 
ble |2] shows one such implementation of the controller in pseudo-C. One of the 
interesting aspects of the processes in this report is their concurrency, which 
can make the structure of the state transitions rather complicated. However, a 
close inspection of the programs reveals a relatively simple transition structure. 
Our investigation thus begins with a forward analysis of both programs. Since 
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/* variable declarations */ 
double Ac [2] [2] ; 
double Cc [2] ; 
double Be [2] ; 
double Dc; 

double xc[2], yc,y,u; 
int flags, flaga; 

/* aux vars */ 
int i,j; 

double xc_new[2] ; 
/* code */ 

Ac [0] [0] =0 . 4990 ; Ac [0] [1] =-0 . 0500 ; Ac [1] [0] =0 . 0100 ; Ac [1] [1] =1 . 0000 ; 

Cc[0]=564.48;Cc[l]=0.0; 

Bc[0]=l. ;Bc[l]=0. ; 

Dc=-1280. ; 

xc[0]=0. ;xc[l]=0. ; 

receive (y) ; 

while(l) { 

yc = (y > 1. ? 1. : y) < -1. ? -1. : y; 

u = 0.0; 

for(i=0;i<2;i++) u += Cc[i]*xc[i]; 
u += Dc*yc; 
for(i=0;i<2;i++) { 
xc_new[i] = 0.0; 

for(j=0; j<2; xc_new[i] += Ac [i] [j] *xc [j] ; 

xc_new[i] += Bc[i]*yc; 

} 

f or (i=0 ; i<2 ; i++) xc[i] = xc_new[i] ; 

send(u) ; 

receive{y} 

} 

Table 2: Implementation of controller in C 
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only local stability has been proven from the program specifications, we may 
not expect any better result by inspecting and documenting the code. Declared 
but non-initialized variables v will be assumed to belong to the empty set, that 
is V G_L. Likewise, the set of all possible values will be denoted T, and a vari- 
able V that may take any possible value will be characterized as w € T. Since 
both programs run concurrently, tracking state values and the sets they belong 
to may be somewhat challenging. The only variable whose value is assumed 
to be initialized prior to the plant execution is the initial plant state Xp. For 
the purpose of stability analysis, we assume that the state Xp is initially within 
the ellipsoid £q, where Q is defined in Thus, prior to the initiation of 

the plant process, and using insight from the system specification analysis, we 
form the condition {(xc^Xp) G £p, {y,yc,u) G T}. The first three lines of the 
plant dynamics Ip, 2p, 3p do not change this condition, and so they may be 
commented as 

{{xc,Xp) e £p, {y,yc,u) g T} 

Ip: Ap = [1. 0000, 0.0100;-0. 0100,1. 0000] ; 

{{xc,Xp) e £p, {y,yc,u) e T} 
2p: Cp=[l,0]; 
{{xc, Xp) e £p, {y, Vc, u) e T} 
3p: Bp = [0.00005; 0.01] ; 

{{xc, Xp) e £p, {y, Vc, u) e T} 

Line4p: while ( 1 ), together with line 9p : end ; defines the main loop of 
the plant dynamics. We use the requirements analysis to postulate the following 
set of pre- and post-conditions 

{{xc, Xp) e £p, {y, Vc, u) e T} 

4p: while (1) 

{{xc, Xp) e £p, {y, yc, u) e T} 

{{xc, Xp) e £p, {y, yc, u) e T} 

9p : end 

{false} 

The last assertion {false} simply indicates a never ending loop; this should 
be expected from a dynamical system that never ceases to execute. In the follow- 
ing, and in order to simplify the notation, we will omit to mention the variables 
whose value could be arbitrary (that is, the variables whose values belong to 
T). With this convention, the previous set of commented lines becomes 
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{{xc,Xp) e £p} 

4p: while (1) 
{{Xc, Xp) e £p} 

{{xc,Xp) e £p} 

9p : end 
{false}. 

The following line, 5p : y = Cp*xp assigns a value to the variable y. The 
corresponding pair of conditions is therefore 

{{xc,Xp) e £p} 

5p: y = Cp*xp 

{{xc,Xp,y) G Gr}, 



The next line, 6p: senci(y) ; does not affect the state variables of interest 
and therefore becomes 

{{Xc,Xp,y) e Gr} 
6p: send(y) 

{{xc,Xp,y) e Gr] 

At this point, the execution of the plant process is blocked by the receive 
command in line 7p. We now turn our attention to the controller code. The 
first four lines of the controller code Ic , 2c , 3c , 4c, are commented using the 
statement {{xc,Xp,y) e Gr, {yc,u) e T} or, in short, {{xc,Xp,y) e Gr}- Thus 
we get 

{{Xc,Xp,y) e Gr} 

Ic: Ac = [0.4990, -0.0500; 0.0100,1.0000]; 

{{Xc,Xp,y) G Gr} 

2c: Cc= [564.48 0] ; 

{{xc,Xp,y) e Gr} 

3c: Bc=[l; 0] ; D = -1280; 

{{xc,Xp,y) € Gr}, 

4c: xc = zeros (2,1); 

{{xc,Xp,y) e Gr}- 

At this point, the execution of the controller program is halted by the com- 
mand 5c: receive (y) and can proceed only after the plant has executed line 
6p. Valid pre- and post-conditions are 



with 




{{xc,Xp,y) G Gr} 
5c: receive (y); 
{{xc,Xp,y) G Gr} 
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We then encounter line 6c: while (1) which, together with line 12c: 
end defines the main loop of the controller. We postulate the following pre- and 
post-conditions for these instructions 

{{xc,Xp,y) e Gr} 
6c: while (1) 

{{Xc,Xp,y) e Gr} 

{{xc,Xp,y) e Gr} 
12c: end 
{false}. 

The next line 7c: yc = max (min (y , 1) ,-1), applies the nonlinear satura- 
tion operator to y. Common systems practice (as well as the methods used 
to compute P in the first place) suggests that this nonlinear relation may be 
accurately captured by the quadratic inequality 

{yc-y){yc-0.2y)<0, 

also named a sector bound on the saturation operator. This sector bound is not 
true in general, but it holds for all variables in Gr- This bound may also be 
expressed in matrix form 





T 


' 













X(^ 




Xc 


T 


Xq 


Xp 


















Xp 




Xp 


T 


Xp 


y 










0.2 


-0.6 




y 




y 


y 


Vc 










-0.6 


1 




yc 




Vc 




yc 



< 0. 



(3) 



Thus one possible set of pre- and post- conditions for line 7c is 

{{xc,Xp,y) e Gr} 
7c: yc = max(min(y, 1) ,-1) 



{xc,Xp,y) e Gr, 





T 




> 


Xq 




Xq 




Xp 


T 


Xp 


< > 


y 




V 










> 



However, following the principles outlined in |FA08] . we seek to capture all 
variables within a single quadratic constraint, for the purpose of making proof 
checking simple and mirroring the usage of the iS-procedure to establish stability 
proofs at the specification level. To do so, we rely on the following lemma: 

Lemma: Assume the real vector variables z and w satisfy the constraints 



Z 

U 



> 



for a given U = U'^ and 





T 




z 


T 


z 




w 




w 



< 



9 



for T = . Then for any real /i, we have 



z 
w 



eGv, V = 





■ / 


■ 




■ u 


■ 


.) 


-1 


■ u 


" 


( 








- M 





/ 







/ 



whenever V , defined as such, exists and is positive definite. 

Proof: The proof is simple and therefore omitted. 

Applying this lemma with U = R, T defined as in ([s]) and fi — 
yields the inequahty 



-A = -6.76 



V 



> 0. 



y 

Vc 

and the pre- and post-conditions 

{{xc,Xp,y) £ Gb} 

7c: yc = max(min(y, 1) ,-1) ; 

{{xc,Xp,y,yc) e Gv} 

At this point, the variable y is not of use anymore and we may therefore 
release it. Thus we replace the post-condition {{xc,Xp,y,yc) G Gv} of line 7c 
by the post condition {{xc, Xp,yc) G Gw}, with 



W 



/ 
0/00 
/ 



V 



/ 
0/00 
/ 



The assignment 8c: u = Cc*xc + Dc*yc is then easily commented as 

{{xc,Xp,yc) e Gw} 

8c: u = Cc*xc + Dc*yc 



{(^ 
with 



Vc) e Gx} 



X 



I 








1 



W 



I 








1 



The next line of the controller loop 9c : xc = Ac*xc + Bc*yc can be treated 
similarly to line 8c, to yield the triple 

{{xc,Xp,u,yc) G Gx} 

9c: xc = Ac*xc + Bc*yc 

{{xc,Xp,u) e Gy} 
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with 



Y = 



Ac Be 
7 
10 



X 



Ac Be 
7 
10 



Note that the variable jjc is not used anymore and has therefore been released. 
Line 10c does not influence the variables and therefore can be commented as 



{{Xc,Xp,U) G Gy} 

10c: send(u) 

{{xc,Xp,u) e Qy} 

The only remaining line in the controller loop is then line 11c : receive (y) , 
which blocks the controller loop until a new value of y is available. 

At that point, it becomes important to look again at the plant process, 
which restarts at line 7p: receive (u). A necessary post-condition for this 
line is {{xc,Xp,u) € Gy}, but no pre-condition can be clearly extracted. Line 
7p will therefore be only partially instrumented, replacing the pre-condition by 

the symbol :. In addition, we indicate between brackets which line execution in 
the controller code (line 10c) unlocks line 7p, yielding the "triple" 



7p: receive (u) [10c] 
{{xc,Xp,u) e Gy}- 

Consider then line 8p : xp = Ap*xp + Bp*u;. This assignment is properly 
documented with pre- and post- conditions to yield 

{{xc,Xp,u) e Gy} 
8p: xp = Ap*xp + Bp*u; 
{{Xc, Xp) e Gz} 

where 

[70 
[0 Ap Bp 

and u has been released. 

At that point, the plant dynamics meets line 9p: end. That line has al- 
ready been instrumented with candidate pre-and post conditions, which there- 
fore must be shown to be compatible with the post-condition of line 9p. More 
precisely, we must show 

{xc, Xp) eGz {xc, Xp) e Ep. (4) 

Numerical computations show that this is indeed true (modulo floating point 
operation arithmetic errors). While the entire loop describing the plant dy- 
namics has been investigated, such is not the case of the controller dynamics, 
whose lines 11c and 12c have not been executed yet. We therefore need to keep 



Y 



I 

A. 
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tracking the execution of the two programs. While the control program is still 
blocked at line 11c, the plant dynamics loop sends the program execution to line 
4p, 5p, and 6p, where the pre- and post- conditions are already established and 
coherent. Line 7p then blocks any further propagation of the plant dynamics, 
while the controller program is now allowed to proceed past line 11c, for which 
a valid post-condition is therefore {{xc,Xp,y) G Gu} and we will instrument as 
follows 

lie: receive (y) [6p] 

{{xc,Xp,y) e Gb} 

which is compatible with the following line, 12c, which has been already 
documented. 

The resulting commented programs then look like those shown in Table |3] 



2 Sector-bounded nonlinearities and stability anal- 
ysis 

Following similar developments in [FAQ_8j , there are apparent differences between 
the stability analysis performed in the requirements-level stability analysis, sum- 
marized by Eq. ([2]), and the forward analysis performed afterwards. Namely, 
we want to show that the invariance condition Q is satisfied if and only if the 
condition Q is satisfied. For that purpose, we begin with the forward analysis 
of the system, beginning at line 7c: yc = max(min(y , 1) ,-1). The assertion 



{Xc,Xp,y) e Gr, 



with 
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°/o Controller Dynamics 

{{Xc,Xp,y) e Gr} 

Ic: Ac = [0.4990, -0.0500; 0.0100,1.0000]; 

{{Xc,Xp,y) G Gr} 

2c: Cc=[564.48, 0] ; 

{{xc,Xp,y) e Gr} 

3c: Bc=[l;0] ; D = -1280; 

{{Xc,Xp,y) e Gr} 

4c: xc = zeros (2,1); 

{{xc,Xp,y) G Gr} 

5c: receive (y) 

{{xc,Xp,y) G Gr} 

6c: while (1) 

{{xc,Xp,y) G ^ij} 

7c: yc = max(min(y, 1) ,-1) 

{{xc,Xp,y,yc) e Gv} 

skip 

{{xc,Xp,yc) e Gw} 

8c: u = Cc*xc + Dc*yc 

{{xc,Xp,u,yc) G Gx} 

9c: xc = Ac*xc + Bc*yc 

{{Xc,Xp,u) e C/y} 

lOc: send(u) 
{{Xc,Xp,u) G Sy} 

11c: receive (y) [6p] 

{{xc,Xp,y) G Sk} 
12c: end 
{false}. 

Table 3: Commented close 



°/o Plant Dynamics 

{{Xc,Xp) G 5p} 

Ip: Ap = [1. 0000, 0.0100; -0.0100,1. 0000] ; 

{{xc,Xp) G £^p} 
2p: Cp=[l,0]; 

{{xc,Xp) G £^p} 

3p: Bp = [0.00005; 0.01] ; 

{(Xc, Xp) G 5p} 

4p: while (1) 

{{Xc,Xp) G fp} 

5p: y = Cp*xp 

{{xc,Xp,y) G Gr} 
6p: send(y); 
{(a;c, a;p, y) G ^ij} 

7p: receive (u); [10c] 

{{Xc,Xp,u) G t/y} 

8p: xp = Ap*xp+ Bp*u; 

{{xc,Xp) G Gz} 
skip 

{(aJcXp) G fp} 
9p : end 

{false} 



[-loop control program 
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Long and somewhat tedious computations indicate that W, where W is 
obtained by erasing the row and column of V corresponding to the the variable 
y, satisfies 



W 



0.6P-iC^ 



p- 



1 





-1 



/ 
0.6C 1 



/ 0.6C^ 
1 



with P — P — O.lGnC^C. Following the forward analysis and substituting 
matrix expressions for their values eventually yield 



Z= A B 



I 
0.6C 1 









1 



I o.ec^ 







1 



[A Bf. 



The invariance condition Q therefore becomes 

p-i 



[A B] 



I 
0.6C 1 







1 







1 



[A P ]^ <p-i 



Using Schur complements a first time yields the inequality 
-P-i [ A B ] 



[30pt] 



I 
0.6C 1 



I 
0.6C 1 



[A B]- 



-P 
n 



< 0. 



Using Schur complements a second time then yields the inequality 

1 T 



I 

0.6C 1 



[A B Y P[ A B ] 



'I ■ 




■ P " 


0.6C 1 




-fi 



< 0. 



Multiplying this inequality to the left and right by 

1 T 



I 

-0.6C 1 



and 



respectively, yields 

[ A B ]^P [ A B ] 



' p 


■ 










+ M 




-A. 
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/ 
0.6C 1 



0.2C^C -0.6C^ 
-0.6C 1 
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3 Conclusions 



This report describes an approach to documenting control programs, whereby 
the control program code is annotated with logical expressions describing the 
set of reachable program states. This approach constitutes the application of 
the Floyd-Hoare paradigm to control programs. It is shown that the the exten- 
sive domain knowledge gathered by control theory about control system speci- 
fications is readily applicable to develop stability and performance proofs of the 
corresponding control programs. Some system elements, such as sector-bounded 
nonlinearities, can be handled via forward analysis, and this analysis was shown 
to be equivalent to well known results in absolute stability theory. 

The analyses discussed in this report may be used in several different con- 
texts: First, they may be used in an autocoding environment, whereby diagram- 
based specification languages such as Simulink can be autocoded in a target 
language such as C to automatically produce software with extensive, inlined 
proofs of software stability and performance. Such a documented software may 
then be easily verified independently by a proof checker. 

Second, the analyses described in this paper may also be used to produce 
stability and poerformance proofs of undocumented software, even if little or 
no information is available about the software specifications. Producing such 
proofs becomes an integral part of the software verification proc;css. Although 
time-consuming, this process is the only available option for much of the existing 
body of real-time control software. A particular case of interest is concerned 
with autocoded sofftware, where prior understanding of the autocoding process 
may help speed up the verification task. 
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